/*  Copyright (C) 2003-2011 Thomas P. Scott and Myron Zimmerman

    Thomas P. Scott <tpscott@alum.mit.edu>
    Myron Zimmerman <MyronZimmerman@alum.mit.edu>

    This file is part of HoloStor.

    HoloStor is free software: you can redistribute it and/or modify
    it under the terms of the GNU General Public License as published by
    the Free Software Foundation, version 3 of the License.

    HoloStor is distributed in the hope that it will be useful,
    but WITHOUT ANY WARRANTY; without even the implied warranty of
    MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
    GNU General Public License for more details.

    You should have received a copy of the GNU General Public License
    along with HoloStor.  If not, see <http://www.gnu.org/licenses/>.

    Parts of HoloStor are protected by US Patent 7,472,334, the use of
    which is granted in accordance to the terms of GPLv3.
*/
/*****************************************************************************

 Module Name:
	arith.c

 Abstract:
	Implementation of arithmetic support routines.

--****************************************************************************/

#include "arith.h"

#ifdef	__KERNEL__
#include <asm/div64.h>
// Code expected by the runtime by code generated by GCC:
// return the quotient of the unsigned division of a and b
unsigned long long
__udivdi3(unsigned long long a, unsigned long long b)
{
	do_div(a, b);
	return a;
}
#endif

#ifdef	_MSC_VER
#define snprintf	_snprintf
#pragma warning(disable:4996)	// C4996: _snprintf: This function may be unsafe
#endif

#define	MAXFMTS		6	// maximum PercentE() calls per printf()
#define	MAXWIDTH	20	// maximum formatted width permitted

//
// Brute force implementation of %e (with precision specification)
// for the ratio x/y.  This is useful for kernel environments where
// floating point and standard libraries are not available.
// Constraints: x >= 0, y > 0, 9 >= precision >= 0.
//
char *PercentE(cnt_t x, cnt_t y, unsigned precision)
{
	static char BufferPool[MAXFMTS][MAXWIDTH];
	static int wiper = 0;
	char *buffer;
	unsigned significand, whole, frac;
	unsigned i, scale;
	int exp;

	if (x < 0 || y < 0)
		return NULL;					// not supported
	// allocate a buffer to use
	buffer = BufferPool[wiper];
	wiper = (wiper + 1) % MAXFMTS;
	// determine the scale up factor based on precision
	scale = 1;
	for (i = 0; i < precision; ++i)
		scale *= 10;
	// determine the exponent needed to scale up
	exp = 0;
	while (x < y*scale) {
		x *= 10; --exp;
	}
	// assert: y*scale <= x
	while (x >= y*scale*10) {
		y *= 10; ++exp;
	}
	// assert: y*scale <= x < y*scale*10
	exp += precision;					// adjust for scale up
	// round 
	significand = (unsigned)( (x + y/2) / y );
	if (significand >= scale*10) {
		significand /= 10;
		exp += 1;
	}
	// assert: scale <= significand < scale*10
	whole = significand/scale;			//        1 <= whole < 10
	frac = (significand - whole*scale);	// scale/10 <= frac  < scale
	if (precision != 0) {
		char sfmtX[] = "%d.%0Xde%+03d";
		sfmtX[5] = '0' + precision;		// replace the 'X'
		snprintf(buffer, MAXWIDTH, sfmtX, whole, frac, exp);
	} else {
		static const char sfmt0[] = "%de%+03d";
		snprintf(buffer, MAXWIDTH, sfmt0, whole, exp);
	}
	return buffer;
}

#ifdef USE_FLOAT
#ifdef	_MSC_VER
double sqrt(double);			// Disable sqrt() as an intrinsic
#pragma function(sqrt)
#endif

static double sqrt(double x)	// Use a cheap local implementation
{
	double result;
#ifdef	__GNUC__
	__asm__ ("fsqrt" : "=t" (result) : "0" (x) );
#endif
#ifdef	_MSC_VER
	_asm fld x
	_asm fsqrt
	_asm fstp result
#endif
	return result;
}
#endif	// USE_FLOAT

void MetricsInit(metric_t *metric)
{
#ifdef USE_FLOAT
	int i;
	for (i = 0; metric[i].text != NULL; ++i) {
		metric[i].sum = metric[i].sumSq = 0;
		metric[i].nSamples = 0;
	}
#endif
}

void MetricsPrint(metric_t *metric)
{
#ifdef USE_FLOAT
	int i;
	for (i = 0; metric[i].text != NULL; ++i) {
		double mean, rms;
		if (metric[i].nSamples == 0)
			continue;
		mean = metric[i].sum / metric[i].nSamples;
		rms = sqrt(metric[i].sumSq/metric[i].nSamples - mean*mean);
		printf("%s : %.2e +- %.2e (%d samples)\n",
			metric[i].text, mean, rms, metric[i].nSamples);
	}
#endif
}

void MetricSample(metric_t *thisMetric, cnt_t x, cnt_t y)
{
#ifdef USE_FLOAT
	double value = (double)x/(double)y;
	thisMetric->nSamples++;
	thisMetric->sum += value;
	thisMetric->sumSq += value*value;
#endif
}

void arith_self_test(void)
{
#ifdef USE_FLOAT
	printf("sqrt(25.0)=%d\n", (int)sqrt(25.0));
#endif
}
